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Abstract. The transmission of forces through a disordered granular system is 
studied by means of a geometrical-topological approach that reduces the granular 
packing into a set of layers. This layered structure constitutes the skeleton 
through which the force chains set up. Given the granular packing, and the 
region where the force is applied, such a skeleton is uniquely defined. Within 
this framework, we write an equation for the transmission of the vertical forces 
that can be solved recursively layer by layer. We find that a special class of 
analytical solutions for this equation are Levi-stable distributions. We discuss 
the link between criticality and fragility and we show how the disordered packing 
naturally induces the formation of force-chains and arches. We point out that 
critical regimes, with power law distributions, are associated with the roughness 
of the topological layers. Whereas, fragility is associated with local changes in 
the force network induced by local granular rearrangements or by changes in the 
applied force. The results are compared with recent experimental observations in 
particulate matter and with computer simulations. 



PACS numbers: 61.43; 45.70; 62.20 

1. Introduction 



In recent years many experiments 



J,J||J, | @, | |, simulations @ 



and theoretical approaches jl3[ 15f |lq, |17| have been devoted to the study of 
propagation of forces through granular aggregates. One of the interesting features 
that rise from these studies is the formation of strong inhomogeneities in the spatial 
distribution of contact forces (formation of force chains), where a small fraction of 
the grains carries most of the total force (see Fig.|I]). These chains typically extend 
on space scales much larger than the grain dimensions and are associated with broad 
stress distributions |3[ [|, ||, ^[ ||, |l(], Moreover, simulations and experiments show 
strong rearrangements of the force-chains under small changes in the compression axes 
suggesting a 'fragile' behavior associated with the force chain structure [jll], [If], jljfl . 
In the present paper we show that the formation of force chains, the broadening of the 
force distribution and the 'fragile' behavior are direct consequences of the topological 
disorder in the granular packing structure. 

In order to study the static stress distribution, we use a scalar model where the 
force balance is considered in the vertical direction only and torques are neglected. 
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Figure 1. Force chains through a disordered granular packing made by a mixture 
of photoelastic disks with two different sizes. The luminosity is proportional to 
the amount of stress carried 13. 



This approach is similar to several theoretical works already proposed in the literature 
|| [ll], 13, Q and in particular to the well known q model j|, |13| . In our model 
however, we explicitly consider a disordered arrangements of grains (instead of placing 
the grains on a regular array and introducing randomness in the transmission term, 
as for the q model). This leads to strong differences and some important, novel 
consequences (see also Appendix A). 



2. Topological characterization of the granular packing 

The complete description of a granular structure in 2 and 3 dimensions requires 
in general the information about the positions, sizes and shapes of all the grains. 
But, to describe the propagation of the weight (vertical forces), one can use only the 
topological information about the neighbors and the contacts among grains combined 
with the knowledge about the direction and intensity of the transmission among 
neighboring grains (see Fig.||). 

To univocally define the neighboring structure we use the 'Navigation Map' p9| , 
which consists in the construction of space-filling cellular structure with generalized 
Vorono'i Q cells where the faces are made by the set of points equidistant to the 
surfaces of the grains. From this construction two grains are first neighbors if they 
share a common face. The neighboring structure is topologically defined by the matrix 
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C: 

{1 if grain i is first neighbor of grain j, 
otherwise, 

which is symmetric with zeros on the diagonal. This matrix describes the neighboring 
structure of the packing. But, by construction, two neighboring grains i,j (with 
= 1) might not be in contact, whereas two grains which are in contact are 
surely neighbors. Since the force propagates only through grains which are in contact, 
it follows that to study the force-propagation one must also take into account the 
topological structure of the inter-grain contacts (Dodds network) This can be 

done by defining the 'kissing' matrix K: 

if grain i is in contact with grain j, 
otherwise. 

(with Kij — Kj.i, and Kij = 0). The number of grains which are in contact ('Kiss') 
with a given grain i is 

N 

Hi = K *,j , (3) 

3=1 

where N is the total number of grains in the system. An upper bound on the number 
of adjacent grains n iy comes from geometrical considerations: it is indeed known that 
no more that 6 equal disks and no more than 12 equal spheres can stay in touch with a 
central one (kissing numbers) j2^, ^3|. On the other hand, when the sizes can fluctuate 
and the grains can assume different shapes, these numbers might change dramatically. 
But, for grain sizes distributed in a finite range, one can always find an upper bound 
on the number of neighbors in contact with a given grain. For instance, when the 
grain sizes are homogeneous and the shapes isotropic, these upper bounds are very 
close to the above values of 6 and 12, in 2 and 3 dimensions respectively. A lower 
bound on n, is given by the mechanical stability condition, as we will see shortly. 



3. Weight propagation 

We study the propagation of the vertical component of the forces (weights) from grain 
to grain. This is a conservative system: the total weight on the bottom of the stack is 
equal to vertical force applied on the top plus the proper weights of the grains in the 
stack (we neglect the presence of walls that we put at the infinite). 

We describe the propagation of the vertical component of the forces inside the 
granular structure in term of following equation: 

N 

Wi=w° + ^2wjQj,i , (4) 

3=1 

with u>i being the vertical force applied at the barycenter of grain "«" . The quantities 
Qjj are the elements of the weight transmission matrix Q; they are associated with 
the kissing matrix through the relation Qj.i = rij^Kj^i, where the quantities rjj^ take 
values between and 1. These quantities are asymmetric: if r/j.i > then rj^j = 
(but not vice versa), giving in this way the orientation of the force network. 
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Figure 2. A granular packing is topologically described by the neighboring 
matrix (lines). The transmission of the weight is through adjacent grains 

and this structure is the Dodds network which is defined by the 'kissing' matrix 
Kij (arrows). The force propagates in one direction and inhomogeneously among 
adjacent neighbours. This can be described by the weight transmission matrix 
Qi,j — Vi,j • 



The symbol indicates the proper weight of grain i with an exception for the 
top layer of the granular stack where one can apply an external force. In this case, 
the symbol w® indicates the proper weight of grain i plus the portion of external force 
acting on this grain. 

Equation Q is a rather general expression: it assumes only that -as it is- the 
vertical force propagates only through grains which are in contact. Equation |] can 
also take into account non-linear effects if one supposes that the quantities Qj^ are 
dependent on the force ittj. However, in Eq.^ the vectorial nature of the forces and 
the tensorial nature of the stresses have been completely neglected. Nevertheless, we 
will show that the scalar description given by EqJJ is very powerful and can take 
into account complex phenomena like the formation of force chains and arches or the 
fragility. 

In Eq.^| the whole information about the granular packing is contained in 
the weight transmission matrix Q which must therefore be consistent with a real 
connection network in a physically-realizable granular packing. This gives some 
constraints on Q. The conservation of the total weight implies that the sum over 
each column of the matrix Q must be equal to 1, i.e. 



The stability under gravity imposes that each grain must lay on at least 3 other 
grains in 3 dimensions or 2 grains in two dimensions, therefore 

JV 



N 




(5) 




(0) 
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with D the dimensionality of the space and 9(x) the step function (6{x) — 1 when 
x > 6{x) = when x < 0). 

The grain i receives the vertical force from a number 

N 

mi = J2o(Qj,i) (?) 

3=1 

of neighbors, and shares its total weight Wi through 

N 

Si = Y,0(Qi tj ) (8) 
i=i 

neighbors. Clearly, m, + Sj = u, is the number of grains adjacent to i and, from the 
stability condition (Eq.^|), it must be: rii > D, giving therefore a lower bound on this 
quantity. 

4. A linear equation 

We now seek for solutions of Eq.||. Let first consider the special case where the elements 
Qi j of the weight transmission matrix depend on the granular packing structure only 
(i.e. the Qij don't depend on Wj). In this case Eq.[I| is a linear equation, and its 
solution can be written straightforwardly. Indeed, an equivalent way to write Eqj^, is: 

N 

J £w j (5 j , i -Q j , i )=w° i , (9) 
i=i 

with Si j the Kroneker's delta. When the Q^i are independent of Wj, Eq.|| is a set of 
N linear equations with N unknowns. If the determinant of I — Q is different from 
zero, this system will admit a unique solution. Physically this condition corresponds 
to impose that in the force network there are no closed oriented rings. In terms of the 
elements of the force transmission matrix this implies that, for any given i, and for 
any arbitrarily chosen set of sites indexes {ji, jk}, it must hold 

Qi,jlQjl,j2Qi2,j3---Qjk,i = ^ ■ 

Note that this condition implies Q x = for any x > N. Indeed it is impossible 
to make a non-intersecting path with a length larger than the number of grains in the 
system. 

Equation ^ can be also written in a vectorial notation 

w(I-Q)=w° (11) 

and, when the determinant of I — Q is different from zero, its solution is 

w = w°(I-Q)- 1 . (12) 

Although Eq.|l^ is formally trivial, its solution cannot be easily pursued 
analytically in the general case and it might become numerically untreatable for 
large samples. Therefore, in the following we develop an alternative approach which 
considers the system as structured into a set of layers and makes possible to solve Eq.|] 
recursively. 
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5. Construction of the layered structure 

Any granular packing can be formally reduced into a stacking of layers where the 
vertical force components are transmitted downward from layer to layer. The elements 
of these layers can be single grains or sets of grains that collectively behave as a single 
element that receives the weight from the layer above and transmits it to the layer 
below. 

In order to construct such a layered system, let first label all the grains at the 
top of the stacking (on which -eventually- the external load is applied) as belonging to 
the layer zero (Co). Let also consider as belonging to Co all the grains which transmit 
some weight to grains already labelled as belonging to Co- Once defined Co, the first 
layer (C±) is made by the set of grains which are neighbors of grains in Co and don't 
belong to Co- Moreover, are belonging to C\ all the grains which are not yet classified 
in Co and that transmit some force to grains already classified in C\. Analogously, 
in general, the layer C a is constituted of grains (i) that are neighbors of grains (j) 
in layer C a -\ and don't belong to any already classified layer (i 6 C a if i £ Cp with 
(3 < a, and if exists at least one j g C a -\ such that Cjj = 1). In addition, to C a 
belong all the non-classified grains which transmit some force to grains which have 
been already classified in C a (i G C a if i ^ Cp with (3 < a, and if exists at least one 
j G C a such that rjij > 0). 

In the system of layers that we have now constructed we have two kinds of grains: 
(i) some grains receive the weight from the layer above only, and transmit it only to 
the layer below; (ii) other grains receive or transmit weight to grains in the same layer. 
We want to reduce the packing to a set of 'elements' where every element k (which 
belongs to some layer, k G C a ) receives the weight from elements (j) in layer above 
(j G C a -i) and transmits the weight to elements j' in the layer below (j' G C a +i). 
To this purpose we identify as single elements of the layered structure the grains of 
type (i), and the clusters of grains of type (ii). 

A schematic graphical representation of the procedure to identify the granular 
layers and their elements is given in Figure [| 



6. Recursive solutions in layered packings 

In the layered system that we have constructed in the previous section the weight is 
transmitted layer by layer from the top of the staking (layer Co, which receives the 
external load) down to the bottom layer Ch (layer at distance h from the top one). 

In analogy with Eq.||, it follows that the weight on a given element in layer C a 
(a > 1) is given in terms of the weights of the elements on the layer above C a -\ 
through the following equation 

Wi = w® + w jQj,i with i G C a . (13) 

Note that in this case the indices and 'j' refer to elements of the packing instead 
of grains. The definition of Qij is identical to the one given before, but it refers to 
a granular structure which has been reduced to a system of elements disposed into 
layers which make an SSI network |24| (see Fig.||c). The forces on the top layer are 
known: W{ = (for i G Co). Therefore, by using recursively Eq.[ll|, one can calculate 
step by step the forces Wi from the first layer to the bottom. This recursive strategy 
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Figure 3. (a ) The layer Co is made by the set of grains on which the force is 
applied plus all the grains that transmit force to them. C a is made of the grains 
which are neighbors to C a ~ 1 plus the grains which transmit some force to these 
grains in C a but have been not yet classified, (b ) Clusters of adjacent grains 
in the same layer which transmit the force among the mse lves are considered as 
single 'elements', (c) The resulting structure is an SSI |24j layered packing made 
of elements which receive some weight from the layer above and transmit it to 
elements in the layer below. 



can be applied to solve Eq.|l3| even in the non-linear case, with the only assumption 
that the quantities must depend only on the forces applied from the layers above. 



7. Force chains, arches and fragility 



In experiments and computer simulations of force propagation through granular 
matter an important fact is observed: only a small part of grains carries the most 
part of the weight. The force propagates through 'chains', and these chains generate 
'arches'. The formation of these chains implies that there are some region of the 
packing where the weight tends to concentrate in, and there are other regions where 
the weight tends to diffuse out. 

In order to understand this mechanism in term of our layered structure, let us 
consider the shape of these layers. In a crystalline-ordered packing (an FCC stacking 
of balls, for instance) the layers are made of sets of elements disposed on flat surfaces 
which are parallel to each other (or straight lines in 2D packings). When disorder 
is introduced, these surfaces start to bend and become rough. The global curvature 
will rest zero, but locally 'valleys' and 'mountains' will make a dimpled landscape. 
These rough surfaces are staked one upon the other and therefore the roughening of 
neighboring layers must be similar and the local curvature is expected to be preserved 
for a few layers (see Fig.^J) 

Any given element of the layered packing receives the weight from the elements 
in the layer above and distributes it among the neighbors in the layer below. Because 
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Figure 4. (a) The layers in disordered granular packings are made of sets of 
elements placed on rough surfaces (wiggled lines in two dimensions) that make 
dimpled landscape rich of 'valleys' and 'mountains'. For geometrical reasons, 
successive layers have similar local curvature and the roughening is persistent, 
up to a certain distance, through the layers, (b) The force is concentrated 
(or 'focalized') when transmitted through 'mountains' whereas it is diffused 
('defocalized') in the 'valley' region. 



of the disorder, some elements have more neighbors in the layer above than others. 
Clearly, if we had no structural correlation among the layers, we expect that the 
elements that receive the weight from more elements will result -in average- with a 
larger weight than the others with fewer neighbors above. This mechanism would be 
amplified layer after layer if this property of abundance/deficiency of neighbors would 
locally propagate through the layers. This is indeed the case in granular packings 
where, for geometrical reasons, the roughening of the layers must be persistent through 
a certain number of layers. The persistence length is strictly related with the average 
chain length and gives the measure of the amplification effect. In particular, we have 
that elements in 'valleys' tend to dissipate the force (defocalizing effect) by propagating 
it to a larger number of neighbors with respect to the number from which they receive 
it. Whereas, elements in the 'mountains' tend to concentrate the weight (focalizing 
effect) by transmitting the force to a smaller number of elements with respect to the 
number from which they receive it (see Fig.|j). The formation of arches in an 'hands- 
on' simple experiment [181 made with photoclastic material between two polarizing 
filters, is shown in Figfi . The result of a computer simulation of weight propagation 
trough an SSI structure is shown in Fig.|| p5| . As one can see, most of the stress 
propagates through force- chains which make arches by self-interconnecting into an 
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Figure 5. Force chains resulting from a computer simulation of weight 
transmission through an SSI granular system. 



intricate network. 



The readjustment of a few grains or a changement in the point or direction of 
the applied force strongly affects the layered structure. Indeed, a local changement 
modifies locally the composition and curvature of the layer. This modification 
propagates through the whole system by changing the roughening and the set of grains 
in the layers. The final effect is a modification in the whole structure of the force-chain 
network. This long-range propagation associated with the focalizing/defocalizing 
effects, is responsible for some of the critical effects in granular matter and it is at the 
origin of the fragility in these systems. 
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8. Force distribution: an analytical approach 

In this section we calculate analytically, in the linear case, the probability distribution 
for the vertical force components. The probability P(w, r) to find a force w at the 
site r is related to the probability that the grains at sites j have weights Wj , and that 
the propagation of these weights through the contact matrix Qj r produces exactly w, 
which formally writes (from Eqj^): 

P(w,r) =TT / dQ k , r Q(Q) / dwjP(wj,j)S(w - w° - ^Wk'Qk'.r), 
k,j J ° J v 

(14) 

where 8{x) is the Dirac delta function and the quantity f2(Q) is the probability of a 
given force transmission matrix Q. 

Let us introduce the Fourier transform of P(w, r) and its inverse 

PicpJ) = N [ dwe iw *P(w,j) 



P(w,j) = j dpe-^P&J) , (15) 

with J\f a constant which takes into account the normalization to 1. By using this 
transformation, Eq. 14 becomes 

P(w,r) = AfJ] J dQ k , r n{Q) J dw 3 J d^e- iw ^P{^j,j) 

kj 

d( p e -i( w - w ° r -j: k ,w k ,Q k , ^ ^ 

where the Dirac delta function has been written as 

S(x) = TV J d^e""* . (17) 
Equation [l6| can be re-written as: 

P(w,r) = N fd0e-« w - w '> 11 f dQ k ,MQ) II / *PiH<Phi) 

dw j e~ iw ^e iw ^ r4> , (18) 

where the integral over the Q's is restricted only to the set T r of elements which are 
adjacent to the element r and are transmitting the force to this element (i.e. Qj, r ^ 0). 
The integration over Wj yields 

P(w,r) =M [ ^e-^-^TT-^ f'dQk.MQ) IT p (Q 

J k " J ° jeXr 

(19) 

from Eqjl^, its Fourier transform is 

P(^,r) =N 2 [ dcp { dwe iw ^e- i{w - w °^ 



i r 1 

Ht? / dQ k>r n(Q) TT P(Qj, r <t>,j) ■ (20) 

k ^ J ° jEXr 
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The integration over w gives finally 

Pty, r) = JVe iw °^ TT T7 / dQk,M<$ J[ P{Q 3 ^, j) . (21) 

k ^ J ° j£l r 

If we assume that f2(Q) factorizes: 

^({Q})=n^(^) > ( 22 ) 

Eq.^l] reduces to 

P^,r)=Afe iw > [] 77 f <KA{Q)P{tii>,i) ■ (23) 

9. Levi-stable force distributions in a layered-structured system 

When the force propagates downward in a layered system, an element in a layer a 
will receive the weight only from a set of elements in the layer above (a — 1) and will 
transmit it to a set of elements in the layer below (a + 1). Let us suppose that the 
probability distribution of the force depends only on the layer number a. This allows 
to substitute into Eq.p3| the quantity P(ip, r) with P a (ip) 



where is the proper weight of the elements in layer a. 

Supposing that the load on a given element in the layer a — 1 is uniformly 
distributed among the n\_ x neighbors in the layer a, the term flj(Q) becomes: 
CljiQ) = S(Cj ~ lA4-i)' This simplifies EqM to 



-P a ^ 



(25) 



where n~ is the number of grains in layer a — 1 from which an clement in layer a 
receives the force. One can directly verify that Eq.^5] is satisfied by choosing P a (ip) 
in the form 

a 

P a (<ip) = exp(i^ — w^V)exp(-^cr Q |V'| 7 + i/3a a tan y — J V#l 7_1 ) , 

(26) 

with fj, and 7 arbitrary numbers and f3 = for 7 = 1 or /? arbitrary for 7^1. The 
other parameters in Eq.|26| are: 

a 

^=n<«-i) -7 . ( 27 ) 
i/=i 

and 

c « = TT— ■ ( 28 ) 

Expression |2TJ might seem complicated, but it is a broadly studied expression: it 
is the Fourier transform of a symmetric Levy-Kintchine stable distribution (26). It 
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describes a probability distribution which is peaked around the value J21=o c \ w< \l c a- 
It has a width proportional to c 1 / 7 , and it is asymmetric with a skewness proportional 
to [3 G [— 1, 1] ([3 positive distribution skewed to the right, (3 negative distribution 
skewed to the left). The parameter 7 G (0,2] is called characteristic exponent. For 
instance, when 7 = 2 the distribution is a Gaussian, whereas for 7 = 1 and (3 = 
it is a Lorentzian. These distributions are stable with respect to the sum of a set of 
stochastic variables and therefore appear naturally in the context of the Central Limit 
Theorem. For 7 < 2 they have tails which decrease as power laws with exponents 
that tend to — (7 + 1) for very large deviations. In general the moment of order n 
becomes infinite when n > 7, it turns out therefore that the variance is only defined 
when 7 = 2, whereas the average is defined for 7 > 1. 

From Eqs.|26|, |27] and |2^ two opposite scenarios come naturally out depending 
on the values of the parameters n + and n~ through the layered structure. Let 
first consider the case where n~ > r£"-i (^ ne elements receive the weight from a 
number of elements which is larger than the number of elements to which the force is 
transmitted). From Eqs. p7| and |2§| it follows that in this case the sum ^a=o c \ w \/ c a 
tends to diverge with a and consequently the value of a tends to infinite. This is the 
mathematical description of the formation of a force chain where, layer after layer the 
force in the chain increases (focalization effect) with an associated broadening in the 
force distribution which becomes a power law with average and spread that tend to 
infinite. The opposite scenario corresponds to the case when n~ < ri^_ 1 (the elements 
receive the weight from a number of elements which is smaller than the number of 
elements to which the force is transmitted). In this case, Eqs. and [28] yield to 
a force-distribution which has a — > and concentrates around the proper weight of 
the elements. The elements receive zero weight from above and diffuse their weights 
outward from this region (defocalizing effect). These are the uncharged grains that 
stay in the regions below the arches. 

10. Conclusions 

We have presented a study of the weight propagation through a granular staking which 
is based on a realistic geometrical-topological description of the packing structure. 
The guiding idea was to analyze the disordered packing as structured in layers 
made of elements which receive the weight from the layer above and transmit it to 
the layer below. This system is an SSI packing and the equation for the vertical 
forces can be solved recursively. We show that the appearance of inhomogeneous 
force distributions into networks of force-chains and arches are naturally found as 
a geometrical consequence of the rough shape of the layers in disordered packings. 
The force is focalized into chains or diffused out and defocalized depending on 
the local curvature of the layered structure. Strong changements in the force- 
network structure can be generated by small rearrangements in the packing which 
are amplified by the geometrical correlations through the layered system. This leads 
to an intrinsic fragile behavior of these systems. We find a class of solutions for the 
distribution of forces which falls in the class of Levy-Kintchine stable distributions. 
The focalizing/defocalizing effects are analytically retrieved for these solutions and 
associated with the local topological properties (n + , n~) of the layered structure. 

Acknowledgments: many thanks to Mario Nicodemi for fruitful discussions. 
T.A. acknowledges many discussions with N. Rivier and C. Oguey. 



Stress transmission in granular matter 



13 



Appendix A. Special cases: the q-model and lattice-based models 

In the literature several models have been presented for the force propagation in 
granular matter. Most of these models assume that the granular packing is stacked 
in an ordered lattice-like structure and the disorder is introduced by supposing an 
inhomogeneous propagation of the force among neighboring grains. This is -for 
instance- the working- framework of the well-known q-model [ll|. It might be 
important to point out that Equation ^ becomes identical to the one of the q-model 
if: 

• the equation is linear {Qij independent on Wj); 

• the granular packing is structured in layers; 

• the force propagates only downward from layer to layer; 

• the constraint given by Eq.^| is satisfied with the equality. 

Therefore, the approach presented in this paper is applicable to the q-model, but 
not vice versa. In particular, in the framework of the q-model, the layer roughening 
(which is the mechanism that leads us to the formation of force chains and arches) 
is associated with the probability that some neighboring grains transmit no force, 
whereas the geometrical correlations between successive layers must be taken into 
account by imposing some 'propagation' of these broken bonds. These are indeed the 
mechanisms, proposed in the literature, which allow, within the framework of a lattice- 
based models (as the q-model), to obtain the formation of arches and critical force 
distributions |l4|]. In our approach these features come directly from the geometrical 
nature of the layer system. 
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